perm filename BBESCO.F4[1,MUS] blob
sn#078092 filedate 1973-12-18 generic text, type T, neo UTF8
00100 C This is taken from the "Handbook of Mathematical Functions"
00200 C by Abramowitz and Stegun, Dover press S1272, 1965.
00300 C Chapter 9 is on Bessel Functions. The polynomial
00400 C approximations are found in section 9.4, pages 369
00500 C and 370. The recursion formula for generating higher
00600 C orders from J0 and J1 is found at the bottom of
00700 C table 9.1, page 390
00800
00900 SUBROUTINE BESCOEF(MI,X,IX,K)
01000 REAL MI,J0,J1
01100 COMMON FREQ(3,0/50,100),FUNC(100),AMP(100),II(1),IJJ(4000)
01200 IF(IX.GT.1)GO TO 30
01300 IF(MI.GE.3.0)GO TO 10
01400 xs=(MI/3)**2
01500 j0=1.0+xs*(-2.2499997+xs*(1.2656208+xs*(-0.3163866+
01600 1xs*(0.0444479+xs*(-0.0039444+xs*0.00021)))))
01700 j1=MI*(.5+xs*(-0.56249985+xs*(0.21093573+
01800 1xs*(-0.03954289+xs*(0.00443319+xs*(-0.00031761+
01900 2xs*0.00001109))))))
02000 GO TO 20
02100 10 xs=3.0/MI
02200 j0=(1.0/sqrt(MI))*(0.79788456+xs*(-0.00000077+
02300 1xs*(-0.0055274+xs*(-0.00009512+xs*(0.00137237+
02400 2xs*(-0.00072805+xs*0.00014476))))))*
02500 3cos(MI-0.78539816+xs*(-0.04166397+
02600 4xs*(-0.00003954+xs*(0.00262573+
02700 5xs*(-0.00054125+xs*(-0.00029333+
02800 6xs*0.00013558))))))
02900 j1=(1.0/sqrt(MI))*(0.79788456+xs*(0.00000156+
03000 1xs*(0.01659667+
03100 1xs*(0.00017105+xs*(-0.00249511+xs*(0.00113653-
03200 2xs*0.00020033))))))*
03300 3cos(MI-2.35619449+xs*(0.12499612+xs*(0.0000565+
03400 4xs*(-0.00637879+xs*(0.00074348+xs*(0.00079824-
03500 5xs*0.00029166))))))
03600 20 FREQ(2,0,K)=J0
03700 FREQ(2,1,K)=J1
03800 FREQ(2,2,K)=J1
03900 RETURN
04000 30 Y=IX
04100 IF(MI.LE.0.0000001)GO TO 50
04200 IJ=IX*2
04300 X=((2.0*(Y-1.0))/MI)*FREQ(2,IJ-2,K)-FREQ(2,IJ-4,K)
04400 40 RETURN
04500 50 X=0
04600 RETURN
04700 END